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Abstract: We revisit the production of leptonic asymmetries in minimal extensions of the 
Standard Model that can explain neutrino masses, involving extra singlets with Majorana 
masses in the GeV scale. We study the quantum kinetic equations both analytically, 
via a perturbative expansion up to third order in the mixing angles, and numerically. 
The analytical solution allows us to identify the relevant CP invariants, and simplifies 
the exploration of the parameter space. We find that sizeable lepton asymmetries are 
compatible with non-degenerate neutrino masses and measurable active-sterile mixings. 
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1. Introduction 

One of the interesting potential implications of (Majorana) neutrino masses is the gener¬ 

ation of a matter-antimatter asymmetry in the Universe. It has been demonstrated that 
the generation of sizeable leptonic asymmetries, leptogenesis, is generic in extensions of the 

Standard Model that can account for neutrino masses [1]. In particular two new ingredients 

are essential for this mechanism to work: the existence of new weakly interacting particles 
that are not in thermal equilibrium sometime before the electroweak phase transition and 
the existence of new sources of CP violation. 

Leptogenesis from the out-of-equilibrium decay of heavy Majorana fermions that ap¬ 
pear in type I seesaw models [1] has been extensively studied (for a comprehensive review 
see e.g. [2]). The simplest version requires however relatively large Majorana masses > 10 8 
GeV [3, 4] (or > 10 1 * * * * 6 if flavour effects are included [5]), which imply that this scenario 


- 1 - 




would be very difficult to test experimentally. It is possible to have sizeable asymmetries 
for smaller masses if a large degeneracy exists, through resonant leptogenesis [6]. 

On the other hand, for Majorana masses in the GeV range, when the neutrino Yukawa 
couplings are small, another mechanism might be at work. In particular, the non-equilibrium 
condition takes place not in the decay, but in the production of the heavy sterile neutrinos. 
The small Yukawa couplings imply that some of the species might never reach thermal 
equilibrium and a lepton asymmetry can be generated at production and seed the baryon 
asymmetry in the Universe. This mechanism was first proposed by Akhmedov, Rubakov 
and Smirnov (ARS) in their pioneering work [7] and pursued, with important refinements 
in refs. [8, 9]. For a recent review and further references see [10]. In most of these works, 
the case of just two extra sterile species is considered, which is also the limiting case of 
the so-called z/MSM where there are three species, but one of them plays the role of warm 
dark matter (WDM) and is almost decoupled, having no impact in the generation of the 
lepton asymmetry. When the mechanism involves just two species, it has been found that 
the observed baryon asymmetry is only possible if the two states are highly degenerate in 
mass. This however was not the conclusion of the ARS paper. 

The purpose of this paper is to explore systematically the parameter space in the case 
of three sterile species (which encompass the one with two neutrinos) as regards the baryon 
asymmetry, in particular we do not want to restrict the parameter space to have a WDM 
candidate. The model has many free parameters (only 5 out of the 18 parameters are 
fixed by the measured light neutrino masses and mixings) and the exploration of the full 
parameter space is challenging. Only with the help of approximate analytical solutions to 
the kinetic equations this task is feasible. The analytical solutions furthermore allows us 
to identify the relevant CP invariants and to reach regions of parameter space where the 
equations become stiff and very difficult to deal with numerically. 

The paper is organised as follows. In section 2 we present the model, which is essen¬ 
tially a generic type I seesaw model, establish the notation and discuss on general grounds 
what are the CP reparametrization and flavour invariants we expect to find in computing 
any CP violating quantity such as any putative lepton asymmetry. In section 3 we present 
the kinetic equations that describe the production of sterile neutrinos and solve them ana¬ 
lytically via a perturbative expansion in the mixing angles up to the third order. In section 
4 we compare the analytical and numerical solutions for several choices of the parame¬ 
ters, and identify the region of parameter space where the analytical solution accurately 
describes the numerical one. In section 5 we use the analytical solutions and perform a 
Monte Carlo scan (using the software package MultiNest [11, 12]) to find regions of param¬ 
eter space that can reproduce the observed baryon asymmetry, and that are compatible 
with the measured neutrino masses and mixings. In section 6 we conclude. 


2. Minimal Model of neutrino masses 

We will concentrate on the arguably simplest model of neutrino masses that includes three 
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right-handed singlets. The Lagrangian is given by: 

3 1 

C = Csm~J 2 L aY ° l ^N' R - Y, 2 ^M lJN n + h -^ 

a,i i,j=l 

where Y is a 3 x 3 complex matrix and M a diagonal real matrix. The spectrum of this 
theory has six massive Majorana neutrinos, and the mixing is described in terms of six 
angles and six CP phases generically. One convenient parametrization for the problem at 
hand is in terms of the eigenvalues of the yukawa and majorana mass matrices together 
with two unitary matrices, V and W. In the basis where the Majorana mass is diagonal, 
M = Diag(Mi, M 2 , M 3 ), the neutrino Yukawa matrix is given by: 

Y = Y f Diag(yi, y 2 ,2/3)kP- (2.1) 

Without loss of generality, using rephasing invariance, we can reduce the unitary matrices 
to the form 1 : 


W = 1/(012,0i3,023, «J) t Diag(l, e ia \e ia2 ), 

V = Diag(l, e^)t/(0i 2 , 0 i 3 ,0 23 , 6), (2.2) 


where 2 


17 ( 01 , 02 , 03 , S) 


( cos 0 i sin 9\ 
— sin 9 1 cos 6 1 
v 0 0 


^ cos 02 0 sin 02 e lS 

0 1 0 

y — sin 02 e 1 * 0 cos 02 


/ 1 0 0 
0 cos 0 3 sin 0 3 
y 0 — sin 03 cos 03 


(2.3) 


Obviously not all the parameters are free, since this model must reproduce the light 
neutrino masses, which approximately implies the seesaw relation: 


-y y M y > 


(2.4) 


where v = 246 GeV is the vev of the Higgs. On the other hand, the known neutrino masses 
and mixings do not give us enough information to determine the Majorana spectrum, not 
even the absolute scale. Very strong constraints can be derived from neutrino oscillation 
experiments for masses below the eV range [13, 14, 15, 16]. Cosmology can exclude a huge 
window below 100 MeV [17, 18, 19, 20, 21, 22], except maybe for one species that could be 
lighter provided the lightest active neutrino mass is below < 3 x 10 _3 eV [20, 21]. The GeV 
range is interesting because an alternative mechanism for lepton asymmetry generation 
could be at work [7, 8 , 9]. Majorana neutrinos in this range are heavy enough to safely 
decay before Big Bang Nucleosynthesis, while they are light enough that they might have 
not reached thermal equilibrium by the time of the electroweak phase transition (EWPT), 
behaving as reservoirs of a putative lepton asymmetry. 

Our goal in this paper is to explore the full parameter space of this model allowed by 
neutrino masses, as regards leptogenesis. An essential condition will be that at least one 

1 Although we use the same notation for the mixing angles and phases of W as those in the usual PMNS 
matrix, they should not be confused. 

2 Note the unconventional ordering of the 2x2 rotation matrices in U. 
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of the sterile neutrinos does not reach thermal equilibrium before the EWPT. This can be 
ensured assuming a large hierarchy in the yukawas [7]: 

2/3 <2/1,2/2- (2.5) 

It is mandatory, however, to have an accurate analytical description, since the uncon¬ 
strained parameter space is huge. We will solve the quantum kinetic equations in a per¬ 
turbative expansion in the mixings in the next section. Since the lepton asymmetry is 
necessarily a CP-odd observable, on general grounds we can derive what are the expecta¬ 
tions in terms of weak-basis CP invariants. 

2.1 CP invariants 

In [23], weak basis (WB) invariants sensitive to the CP violating phases which appear in 
leptogenesis, within the type I seesaw model, were derived. All of them should vanish if 
CP is conserved, and conversely the non-vanishing of any of these invariants signals CP 
violation. They must be invariant under the basis transformations: 

4 W l £ l , 

Nr WrNr . (2.6) 

Defining h = Y^Y, and Hm = M, a subset of the invariants can be written as: 

h = ImTr \hH M M* h*M\, (2.7) 

I 2 = ImTr [hH^M* h* M], (2.8) 

h = ImTr [hH 2 M M*h*MH M \- (2.9) 

Since the are WB invariants, we can evaluate them in any basis. In the WB where the 

sterile neutrino mass matrix M is real and diagonal, one obtains: 

h = MiM 2 AM 2 2 1 Im(/i? 2 ) + MiM 3 AM 3 2 1 Im(^ 3 ) + M 2 M 3 AM 3 2 2 Im(/4), 

1 2 = MiM 2 (M 2 4 - Mf)lm(h 2 12 ) + MiM 3 (M| - M?)lm(h 2 13 ) 

+ M 2 M 3 (M$ - M 2 )Im(h 23 ), 

1 3 = M'i M 2 AM 2l Im(/i 2 2 ) + M 1 3 M 3 3 AM 3 2 1 Im(/i? 3 ) + M 2 3 M|AM 3 2 2 Im(/t| 3 ), 
where AM 4 2 = M 2 — M 2 and, using the parametrization of eq. ( 2 . 1 ) 

Im(h 2 ) = Im[(yty)2 ] = ^ y 2 y 2 Im [W^W^W aj W N ]. (2.13) 

a,/3 

It is explicit in the above expression that such unflavoured invariants depend only on the 
CP phases of the sterile neutrino sector, which are encoded in the unitary matrix W: one 
Dirac-type phase, 5 and two Majorana-type phases au, a 2 . Not surprisingly, these invariants 
are the relevant ones in unflavoured leptogenesis, i.e., in the conventional computation of 
the CP asymmetry generated by heavy Majorana neutrino decay neglecting flavour effects. 


( 2 . 10 ) 

( 2 . 11 ) 

( 2 . 12 ) 
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The combinations of W matrix elements which appear in Im(/i?-) can be expressed in 
terms of the rephasing invariants defined in [24] as follows: 


Im [W^W^W aj W Pj ] 


Im [Wgj Wpj W*j Wfij (Wgj W*j ) 2 ] 

I W ai w aj \ 2 


(2.14) 


Notice that J\y = ± Im [ W ai Wg t W* 3 Wpj} is the Jarlskog invariant for the matrix W, 
while the quantities Im[(WajWY) 2 ] determine the Majorana phases, 0 : 1 , 2 - When consid¬ 
ering processes, such as heavy neutrino oscillations, where the Majorana nature does not 
play a role, only the Dirac phase 6 will be relevant and therefore we expect to find just the 
Jarlskog invariant of the matrix W. 

Since there are six independent CP-violating phases, it is possible to construct three 
more independent WB invariants, which would complete the description of CP violation 
in the leptonic sector. One simple choice are those invariants obtained from /,; under the 
change of the matrix h by h = Y^h(Y, with hi = A^Aj, being \i the charged lepton Yukawa 
couplings, i.e., 

h = Iml \\Y^h e YH M M*Y T h* ( Y*M} , (2.15) 

and analogously for I 2 , J 3 . The corresponding CP odd invariants are Im(/&?■), which in the 
basis where also the charged lepton Yukawa matrix is real and diagonal can be written as: 

Im (hi) = E X l X l Im [Y^Y aj Y Pj Y^\ . (2.16) 

a,(3 


The lepton number (L) violating part of the flavoured CP asymmetries in leptogenesis 
depends on the above combinations [25]: 

4 = E Im [YZiYcjYpjYfcmMi, Mj) , (2.17) 

0,3 

where / is an arbitrary function. Upon substitution of the neutrino Yukawa couplings as 
given in eq. ( 2 . 1 ) can be written as: 

4 = Im [WSiVfaVfJTsjW^WrjMMi, Mj) . (2.18) 

3 PA & 

These asymmetries contain the additional rephasing invariants of the form Im \Wp^p Q Y^ Q Wsj\, 
which depend on the phases in the matrix V(5, (f> 1 , ^ 2 ), showing that the flavoured CP asym¬ 
metries of leptogenesis are also sensitive to the CP phases in the V leptonic mixing matrix, 
besides those in W. 

Alternatively, we choose to construct the WB invariants which will appear when the 
Majorana character of the sterile neutrinos is not relevant, i.e., L-conserving ones. These 
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are given by: 


I[ = Iml \[hH 2 M hH M ] 

= M 1 2 M 2 2 AM 2 2 1 Im(^i 2 /i2i) + M 2 M|AM| 1 Im(/ii 3 ^3i) 

+ Af|A/ 2 AIm( h-i's /i. 32 ) , (2.19) 

J' = 

= M x 2 M 2 2 (M 2 4 - M 1 4 )Im(/ii2^2i) + M 2 M 2 (M 4 - M 4 )Im(/ii 3 /i,3i) 

+ Ml Ml (Ml - Ml)lm(h 2:i h: V2 ) , (2.20) 

J' = 

= M 4 M|AM| 1 Im(/ii 2 /i2i) + Mi 4 M|AM 3 2 1 Im(/n 3 ^ 3 i) 

+ M 4 A'/ 4 AAf. 2 2 Im(/i23/i32) , (2.21) 

where 

MM*) = E . (2.22) 

a, 13 

The L-conserving CP asymmetry in leptogenesis via heavy neutrino decay, as well as 
the CP asymmetries encountered in leptogenesis through sterile neutrino oscillations, are 
sensitive to the above combinations of Yukawa couplings [25]: 

4a = E I m\Y < »K j Yf> j Y$ i \ f(Mi, Mj) , (2.23) 

j,P 

where / is an arbitrary function, and can be written in terms of the rephasing invariants 
as: 

4a = - E E ywyl im/( m*, m 3 ) . ( 2 . 24 ) 

3 /S.^o- 

Notice that the crucial difference between the L-violating and the L-conserving CP asym¬ 
metries is that in ef a the combination of W matrix elements is such that all dependence 
on the Majorana phases a\- 2 disappears, as expected. 

In the approximation of neglecting 7/3 <C yi,y 2 , we obtain that Im \Y ai Y*-(Y'Y)ij\ = 
^Im [YajYZjYpjYfa] reduces to 

Im[Y ai YT(Yty)..] = y 2 y 2 (|F 2a | 2 - | V^ln^W^W^W^il 

+ VW 2 { [yl\W 2i \ 2 - yl\Wu\ 2 ] MWl 3 V la V 2 * a W 2j \ 

+ [yflWijl 2 -ji\W 2j \ 2 ]hn\W^ la V^W 2i ]} , (2.25) 

so in principle we expect that the lepton asymmetry will depend on ten CP invariants, 
namely ImfW^ Vi a V 2 * a W 2 i], with i = 1,2,3 and a = 1, 2, 3 and Jw- 

However, they are not all independent. In ref. [24] it has been shown that in the 
minimal seesaw there are only six independent CP invariants that can be made out of the 
matrices V,W. Two of them correspond to the Majorana phases of W, 01 , 2 , which as 
we have argued before will not contribute in the limit of small sterile neutrino Majorana 
masses that we are considering. Other two are the equivalent of the Jarlskog invariants 
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for the matrices V, W and therefore determine the Dirac phases, 6,6 , respectively. The 
last two are of the form Im [ Wj*Vj a VfaW 2l ] ■ for two reference values of i,a, that fix the 
additional phases fa )2 . Moreover, it can be shown that since we are neglecting the Yukawa 
coupling 2 / 3 , the phase fa of the matrix V does not appear in eq. (2.25), thus we are left 
with only three independent invariants. 

The unitarity of the mixing matrices V, W implies that 

E y ^ = 0, (2.26) 

Q 

E W^W 2i = 0 , (2.27) 


which allows to write the invariants for a = 2 in terms of those with 

a = 1, 3, and the invariants for i = 2 in terms of the corresponding ones with * = 1,3. By 
exploiting the identities 


iin[w; l v l eVfaw 2j \ 


Im[( WWaV&WnHVZpVtoV&Vw)] 
\V la V 2a \i 


(2.28) 


we can write for instance one of the invariants with = 3 in terms of the invariant with 
a = 1 and the Jarlskog invariant for V, Im[ Y 2 a V\ Q V\p\ = ±Jy. 

It is simpler, though, to write the results in terms of the following four invariants, even 
if only three are independent, expanded up to 3rd order in the small mixing angles 9ij, 9ij\ 

/| 2) = -Im^YnV^^] ^ 012012 8^!, 

l[ 3) = Im [W* 2 V 13 V 23 W 22 ] - 012013023 Sin(0 + fa), 

4 3) = Im[W? 3 V 12 V 2 * 2 W 23 ] ~ 0 12 0 13 023 sin((5 - fa), 

J w = -Im[W 2 * 3 W 22 W 3 * 2 W 33 ] ^ 0 12 9 13 e 23 sin<5. (2.29) 


A generic expectation for the CP-asymmetry relevant for leptogenesis is 

A cp = E \ Yak \ 2 A a , (2.30) 

a,k 

with 

A« = E 4 = E MYaiY^YUfiMt, Mj). (2.31) 

* i,3 

Since the CP rephasing invariants are at least second order in the angles, we just need to 
take the diagonal elements in A cp, to keep the result up to 3rd order. Then, in the limit 
2/3 = 0, we get: 


A cp = Vh£(l£ - Vf) E YniW^WfaV^WfafiAfa Mj) 

+ 2 / 12/2 ((2/ 2 2 - 2/i) {/i (2) b(Mi) - g(M 2 )] + I 2 3 \g(M\) - <?(M 3 )]} (2.32) 

- yll i 2\g{M l )-g{M 2 )}), 
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where 


g{Mi) = ylifiM^Mi) - f(Mi , Mi)] - y 2 2 [/(M 2 , M t ) - f{M i} M 2 )\ . (2.33) 


From the above definition of g(Mi), it immediately follows that g{M\) — g(M 2 ) = 
(y 2 “ Vi)[f( M 1 , M 2 ) - /(M 2 , Mi)], so A C p simplifies to 

A C p = ylvlivl - yl) Y i m [WiiWijW£jW 2 i\f(Mi, m 3 ) 



+ ymivl - yl) {\(vl - vl)4 2) ~ yl 1 { 3) 1 [/(Mi, m 2 ) - /(m 2 , m x )] 


(2.34) 


We will see in the next section that this is precisely the yukawa and mixing angle 
dependence we will find when solving the kinetic equations, which is a strong crosscheck 
of the result. 

3. Perturbative Solution of the Raffelt-Sigl equation 
3.1 Sterile neutrino production 

Our starting point is the Raffelt-Sigl formulation [26] of the kinetic equations that describe 
the production of sterile neutrinos in the early Universe. The density matrix is the expec¬ 
tation value of the one-particle number operator for momentum k: pAr(fc) for neutrinos, 
and pN(k ) for antineutrinos. We will assume that only sterile neutrinos and the lepton 
doublets are out of chemical equilibrium, but assume that all the particles are in kinetic 
equilibrium, using Maxwell-Boltzmann statistics: 


Pa{k ) — A a Pe q(fc)j A a — , p a (k) — A a p ec ^{k) , A a — G 


(3.1) 


where p eq (k) = e fc °/ T , with /cq = |k|, and p a denotes the chemical potential normalised by 


the temperature. We will furthermore neglect spectator processes and the washout induced 


by the asymmetries in all the fields other than the sterile neutrinos and lepton doublets. 


We expect this approximation to give uncertainties of 0(1) which for our purpose is good 
enough [27]. 

In [7], only the asymmetry in the sterile sector was considered, neglecting the feedback 
of the leptonic chemical potentials. In this case, the equations get the standard form 


Pn = ~i[H, p N \ - -{r,^ - Peq} 


(3.2) 


and the analogous for pm with H —y H*, where H is the Hamiltonian (we neglect matter 
potentials for the time being but we will include them later on) 



(3.3) 
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r is the rate of production/annihilation of sterile neutrinos in the plasma, which is diagonal 
in the basis that diagonalises the neutrino Yukawa’s: 

r = Diag(ri,r 2 ,0), F ?: oc yf, (3.4) 

where we assume y 3 = 0. In deriving eq.(3.2) it is assumed that the particles involved 
in the production/annihilation of the sterile neutrinos are in full equilibrium (all chemical 
potentials vanish), and that kinematical effects of neutrino masses are negligible. 

Note that only the matrix W appears in these equations and therefore any CP asym¬ 
metry generated can only be proportional to the invariant Jw which depends at third order 
on the mixing angles of W. 

In [ 8 ] it was correctly pointed out that the asymmetries in the sterile sector will be 
modified by the leptonic chemical potentials that will be generated as soon as sterile neu¬ 
trinos start to be produced. Including the evolution of the leptonic chemical potentials has 
two important consequences: new sources of CP violation become relevant and washout 
effects are effective. Leptons are fastly interacting through electroweak interactions in the 
plasma and therefore it is a good approximation to assume they are in kinetic equilibrium. 

An important question is what is the flavour structure of these chemical potentials. 
For T < 10 9 GeV the Yukawa interactions of the tau and muon are very fast, which implies 
that p will be diagonal in the basis that diagonalises the charged lepton Yukawa matrix, 
since no other interaction changing flavour is in equilibrium before the heavy neutrinos 
are produced. Note however that this is not the basis where the neutrino Yukawas are 
diagonal, the two are related by the mixing matrix V. As a result, when the evolution 
of the lepton chemical potentials is taken into account, the CP phases of the matrix V 
become relevant. 

Adapting the derivation of [26] to this situation, we find that the evolution of the 
CP-even and CP-odd parts of the neutrino densities: p± = Pn ± Pn and the lepton chemical 
potentials , p a , to linear order in p, a ,p_, satisfy in this case: 

P+ = —i[H Te , p + \ + [H im , p_] - l^±jk {Y ^Y, P+ - p eq } 

+i'y b N lm[Y ] pY]p eq + {im [Y*pY],p+}, 

P- = -i[H Te ,P-\ + [iW+] - 7 ^ + 7 ^ {yty ;/Q _} 

+^ b N Re[Y^Y] Peq + ^{R ,e[Y^Y),p + }, 

Pa = ~Pa ( 7 ^Tr[YYt /a ] + 7 “ T r[ Re [Yt/ a Y]r + ]) 

+(7“ + it) (Tr[Re[Yt/ Q y] r _] + iTr[lm[Y^ I a Y]r + ]) , (3.5) 

where H re = Re[iL], H ]m = Im[7f], I a is the projector on flavour a and are the 

rates of production/annihilation of a sterile neutrino or a lepton doublet neglecting all 
masses, after factorizing the flavour structure in the Yukawas, 

7 N(l) = YY [ Pe q (Pi)\M { “^lf{2ir)*5(k +P1 -P2 - Ps), (3.6) 

i "Pl 1P2 ,P 3 
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a 


b 




Figure 1: a, b topologies for annihilation/production of sterile neutrinos 


where k is the momentum of the N or v and a(b) refer to the s-channel (t,u-channels) 
depicted in figure 1. In topology a the lepton and sterile neutrino are both in the initial or 
final state, while topology b corresponds to those diagrams where one is in the initial and 
other in the final state. Finally 

Ei f Pl , p 2 ,p 3 p±(pi)\^ifi 2 (^) 4 Hk + Pi -VI- Ps) 

Ei/p 1 , P2 ,p3^q(pi)|A^2 ) | 2 (2vr) 4 (5(fc+pi -p 2 -p 3 ) 

A similar derivation can be found in [28] and we agree with their findings. 

These equations reduce to those in eq. (3.2) in the limit p —> 0 with T ? ; = yf( 7^ + 7^)- 
Most previous studies have assumed that the rates are dominated by the top quark 
scatterings. In this case, the rates are given (in the Boltzman approximation) by the 
well-known result [29, 30] 

1n,q = 2 1n,q = 2 7 1,q = 4 7“,q = ( 3 ‘ 8 ) 

The factor of 2 difference between the rates of the N and the v is due to the fact that the 
lepton is a doublet and the sterile neutrino is a singlet. Note that there is a non-linear 
term of the form 0(pp+), as first noted in [28]. More recently in [31], the equations have 
been written in terms of the PB-L a /3 chemical potentials, however not all the chemical 
potentials (e.g. higgs and top quark) have been included. A full treatment including all 
chemical potentials will be postponed for a future work, but we expect that including these 
spectator effects will change the results by factors of 0 (1). 

In [30, 32], it has been pointed out that the scattering processes LN <->• WH get a 
strong enhancement from hard thermal loops and are actually the dominant scatterings. 
The results of [30, 32] however do not include the chemical potentials of spectators, so it 
is not clear how to include them consistently in the above equations. We will neglect these 
effects in the following. Note however that the lepton flavour structure of these and of the 
top quark scatterings is the same. 
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It is easy to see also that total lepton number is conserved as it should: 

2^/i Q + TY[r_] =0. (3.9) 

Q 

Two approximations are often used in solving these equations: 1) assume that the mo¬ 
mentum dependence of p± follows that of p eq , i.e. kinetic equilibrium for the sterile states, 
which implies r± = p±/p e q are constants and the integro-differential equations become just 
differential equations, 2 ) neglect the kft dependence of the rates by approximating 

(K 1 ) - ^2~- ( 3 - 10 ) 

The effect of these approximations has been studied numerically in [28] and the results do 
not differ too much. We will therefore adopt both approximations that simplify consider¬ 
ably the perturbative treatment. 

3.2 Lepton asymmetries in the sterile sector 

We are going to solve these equations perturbing in the mixing angles up to third order. 
We first consider the simpler case, neglecting leptonic chemical potentials and considering 
in turn the evolution in a static Universe and in the expanding case. 


3.2.1 Static Universe 

We start with eq. (3.2) and assume y% = 0. In this case, neither H nor T depend on time. 
Defining pNij/Pe q = o,%j + ibij and taking into account the hermiticity of p^ we change the 
matrix equation into a vector equation: 


r — (an, « 22 , ai 2 , 612 , 013 , 613 , 023 , 623 , 033 ). 

At 0-th order the system of equations of eq. (3.2) can be rewritten as 

r (0) = A 0 r(°) + h 0 , 

with 

ho = (Tl^eq) r 2Peq j 0, —0), 
and the matrix Aq is constant and has a block structure: 


( (A /) 4x4 0 


An — 


0 


0 (A//) 4x 4 0 | , 

0 0 0 


(3.11) 


(3.12) 


(3.13) 


(3.14) 


Aj = 


/-lb 0 0 

o -r 2 o 


0 0 —Ai 


V 


ri+r 2 

2 


, An = 


/-Ti/2-Ai 3 0 0 

Ais -ri/2 0 0 

0 0 

\0 0 -(A 12 - Ais) - l -f 


Al 2 — A 13 

la 
2 


. (3.15) 
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The matrix can be easily diagonalised and exponentiated so the general solution to the 
equation is 

r ( 0 ) (i) = e Aot f dx e~ A ° x h 0 . (3.16) 

Jo 


At the next order we have to keep 0(8ij ) in the Hamiltonian and translate the matrix 
form into the vector form: 


-i[H^,p^(t)} ->■ Air (0) . 

The equation for the first order correction to the density is 

= AqA 1 ^ + AiA°\t). 


The solution at this order is therefore 

r ( 1 ) (f) = e Aot f dxe- Aox A ir ( 0 Xx). 

Jo 

We can iterate this procedure to get the correction at order n: 

n—l 


rM (*) = A 0 r^ ( t ) + (i), 


i= 1 


with solution 


n—l 


r (n) (t) = e Aot / dxe~ AoX 'VA i A n - i) (x). 

Jo ~1 


i =1 

We can dehne the evolution operator 

U 0 (t,x) = e Aot e~ Aox , 

so that the solution can be written as 

r* 

io 


n—l 


A n \t)= / dx U 0 {t,x)^2A i A n ~ i) (x). 


i —1 


(3.17) 

(3.18) 

(3.19) 

(3.20) 

(3.21) 

(3.22) 

(3.23) 


As a first estimate of the leptonic asymmetry that can be generated, we are interested in 
A/J 33 since this is the sector that will never reach equilibrium (in the absence of mixing) 
and therefore can act as reservoir of the leptonic asymmetry until the electroweak phase 
transition [7]. 

One can easily compute the solution of the eq. (3.21) up to order n = 3, which is the 
first order that gives a non-vanishing result, as expected from general considerations on 
CP invariants. The result at finite t is not particularly illuminating but the limit t —> 00 
is rather simple: 


lim 

t—> OO 


A/?33 

Peq 


t PN33 PN33 

hm - 

t^-oo p e q 


2J\y 


(Ti - r 2 )Ai 2 Ai 3 (Ai2 — A13) 



a 13 ) 2 + 


r 2 " 

1 2 


(3.24) 
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Figure 2: Comparison of numerical (blue) and perturbative (red) solution for A p 33 as a function 
of time, in the case with no expansion of the Universe. The two curves are indistinguishable (left 
plot) until large times (right plot): the vertical line on the lower plot corresponds to (fljjjTi) -1 , 
while those on the upper one correspond to T^" 1 and T)” 1 respectively. 


A few comments are in order. We have not assumed any expansion in F,; in this expression, 
only in the mixing angles. According to general theorems the equations should reach a 
stationary solution if all the eigenvalues of the matrix Aq + A\ + A 2 + ... are real and 
negative. However, because T 3 = 0, one of the eigenvalues of Aq vanishes and it is lifted 
only at second order in perturbation theory, ~ therefore we expect the perturbative 

expansion should break down for t ~ arw > which is the time scale of equilibration of the 
third state. On the other hand, if 8 is small, the perturbative solution should be accurate 
for times t > r~^. Indeed this is precisely what we find comparing the perturbative and 
numerical solutions in figure 2. 

The result is proportional to Jw which is the only CP rephasing invariant that can 
appear in this case. The result vanishes if any two of the masses or the yukawa’s are 
degenerate, since the CP phase would be unphysical in this case. 


3.2.2 Expanding Universe 

Let us turn now to the realistic case of an expanding Universe. As usual, we will consider 
the evolution as a function of the scale factor x = a. in such a way that the Raffelt-Sigl 
equation becomes 


d_ 

dt 


d 

xHu ( x )-Q^p(x,y) 


fixed y 


~i[H{x,y),p(x,y)] - ~{T(x), p(x,y) - p eq {y)}, 


(3.25) 


where H u (x ) is the Hubble parameter, H u = y i?r 45 ^ m?i k ’ anc ^ V = r - Assuming for 
simplicity a radiation dominated Universe with constant number of degrees of freedom, 
during the sterile evolution time we can assume xT= constant that we can fix to be one. 
Therefore the scaling of the different terms is 

A Ad 2 r 1 

H(x,y) = xW — —W t , TiOr) = xH u (x) = —— , (3.26) 

£y x lvi pX 
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where Mp = Mpi anc k y 4 ?r 3 g b (r 0 ) an d 9 *(fao) is the number of relativistic degrees of freedom 
in the plasma during the sterile evolution. 

Therefore the equation as function of x is: 

P = -ix 2 [WAW*,p] - *{7 ,p- Peg}, (3.27) 

where we have defined 

A Ml 

A ij = 7i = (3.28) 

The perturbative expansion works as in section 3.2.1, but now all the A n (x) are x-dependent: 
A n (x) with n > 1 scale like the Hamiltonian, ie. x 2 , while Aq{x) contains terms that scale 
with x 2 and others that do not depend on x. Fortunately, there is an important sim¬ 
plification in that Aq{x) can be diagonalised by an x-independent matrix, therefore the 
path-ordered exponential can be easily evaluated. The result can be written in the same 
form of eq. (3.23), with the evolution operator given by 

Uo(t,r) = e fo A o(o=)dx e -^Aoiy)dy_ (3.29) 


At third order in the mixings, after algebraic simplifications and partial integrations, 
the result can be given in terms of integrals of the form 


Jn{pt It fill Otm finit') — I dx\C 

J 0 


ia 


fXi x 3 

Ja 2 ^+/32X2 


dx 2 e* 


L 


Xn-l x 3 

dx v e ian ^ +pnX '\ 


(3.30) 


where a.i are combinations of Ay and Pi are combinations of 7 j. Up to third order in the 
perturbative expansion only integrals with n < 3 appear. 

Since we are in a regime where 7 j <C |Ay the integrands are highly oscillatory 
and hard to deal with numerically. To evaluate the integrals, we separate the integration 
interval [ 0 ,f] = [0, to] + [to fa] with fa such that folAyj 1 / 3 1 and fo 7 i ^ 1 : 

<Ai (07 j fa j • • j CX-rii Pm t) — Jnfat 1 j fa j •• - Pm ^0) T AJ n (tti , fa , .. -Oini Pm tofa]- (3.31) 

To solve the integrals up to to we can safely Taylor expand in Pi (which results in an 
expansion in 712 /I A*j |) and write the integrals in terms of simpler integrals of the form: 

dxnx'fa e 3 , 

(3.32) 

up to third order in the Pi expansion we just need integrals with n + Yli ki < 3. We can 
use the relation 

4- i F n(fa] = (3-33) 


f k -7LU f xi 
Jnkfal, ..,a n ,t) = / dxixfa e l 3 

Jo Jo 


dx 2 X 3 


k 2 


*2 x 2 

3 


■f 
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with 


/ \ n — z 

F n (a,x) = — 3 s (—ia) 


1+n 

3 r 


1 + n 
3 


1 • 3 

--tax 
o 


(3.34) 


to evaluate immediately the one-dimensional integrals in terms of incomplete T functions. 
The integrals in the range [io,t] can be approximated by the large t behaviour of the 
Ji n (a,t ) functions, after resumming the Taylor series in /3j. Further details are presented 
in appendix A. 

The finite t dependence of the asymmetry A /933 is rather complicated, but the asymp¬ 
totic value is non-zero and rather simple: 


lim 

t—¥ OO 


A /933 
Pe q 


— — Jwll r )2{l2 ~ 7l) 1™ I m [^ 3 o(Ai 2 — A 13 , — A 12 , A 13 , t) + </ 3 o(Ai 2 — A 13 , A 13 , — A 12 , t) 

t —^OO 

+ >^3o(Ai3, —A 12 , A 12 — A 13 , t) + J3o(Ai3, A 12 — A 13 , —A 12 , t)]. 


(3.35) 


This can be simplified to 

j ( A12 Ai 3 \ / A12 A23 \ 

(3.36) 


A/ 9 33 7172(72-71) T 

Inn - = - Jwt a -X X—rr7X Im 

H-OO /9 e q (A13A12A23) 1 / 3 


where 

/ a 9 Aci\ , r°° ■ a X x 3 

^(a^’A^ = ( AiA 2 A3 ) /3 y dxel 3 Jio(A 2j s)Jio(A 3 ,s). (3.37) 

Comparing eq. (3.36) and eq. (3.24) we see that in the expanding case the asymmetry 
is cubic in 7 * and not linear. Note that the dependence on the yukawa’s is precisely that 
expected from a flavour invariant CP asymmetry. In fact this is effectively the situation 
in the expanding case, because the asymmetry is generated at times t <C 7 ~* and the 
dependence in the yukawa’s in this regime is therefore perturbative. This is in contrast 
with the non-expanding case, where the asymmetry evolves all the way till t ~ 7 j 1 . To 
understand the reason behind this different behavior, it is useful to recall the definition of 
A ij from eq. (3.28). Then, we see that Ay-x 3 1 implies AM^/(4T) 3> T 2 /Mp = H U (T), 
therefore in this regime the sterile neutrino oscillations are much faster than the Hubble 
parameter and no asymmetry is produced anymore, since oscillations are averaged out. 
Thus in the expanding Universe the generation of the asymmetry occurs at x ~ |Ajj | -1 / 3 <C 
7f 1 - 

Until now we have neglected the matter potentials, however given the suppression in 
three powers of 7 of the leading result, there are corrections of same order coming from 
the potentials, and in fact they are numerically more important. 

The equation including the potentials in the basis with diagonal neutrino Yukawas is: 

p = -ix 2 [WAW^, p\ -i[v,p] - \{l,p~ Peg}, (3.38) 
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where 


Vij = ^M%Sa = ViS, 


P u ij — 


( 3 . 39 ) 


The result for the asymmetry including the potentials is given by: 

A033 

lim -— = Jw lim^oo Re [z\ J3o(Ai2 — A13, — A12, A13, t) + z 2 </3 o(Ai 2 — A13, A13, —A12, t) 

*->• °° Pe q 

+ ^2«^3o(Ai3, A12 — A13, —A\2,t) + Z3 J3o(Ai 3 , — A12, A i2 — A13, t)] . 

( 3 . 40 ) 


with 


zi = 71 72 A t , + 7iu 2 A 7 + i 
71 V2 - 72^1 + i 


. / 7i7 2 A 0 




22 = 


( 7 i 72 , „ 

Y~ + 2v l y 2 


Z3 = - 7 i 72 A„ - 7 2 uiA 7 + i 
and A v = V2 — v\ and A 7 = (7 2 — 71). 


, ( 7i72A 7 


- 27 iu 2 A, 

A 7 , 

- 27 2 ui At 

-1 


The leading terms 0 (v 2 7) at asymptotic times t > 7 X ) are: 


lim 


A/533 




;JW 


yiym - yt) M l 


t -+00 p eq 2048 vr 3 | AM 2 2 AMf 3 AM| 3 1 x / 3 


( 3 . 41 ) 


( 3 . 42 ) 


where 


« = |Ai 2 Ai3A 2 3| 1 / 3 Im[j3o(Ai 2 — A13, —Ai 2 , A13, t) - J3o(Ai 2 - A13, A13, —Ai 2 , t) 

— <^3o(Ai3 , Ai 2 — A13, — A12, t) + J3o(Ai3, —Ai 2 , Aj 2 — A13, t)] ( 3 . 43 ) 


depends only on the ratios of mass differences and/or the ordering of the states. This result 
is parametrically the same as the result of [ 7 ] if we neglect the dependence of k on the mass 
differences and has the dependence on the yukawas expected from eq. ( 2 . 34 ). 

Considering the naive seesaw scaling yf ~ 2 — , for m v ~ 1 eV and assuming no big 
hierarchies or degeneracies, i.e. Mf AM 2 . M 2 , leads to 


lim 

t^OO 


A/533 

Pe q 



( 3 . 44 ) 


The asymmetry is highly sensitive to the light neutrino mass. Note that we have pushed the 
value to the limit, a light neutrino mass in the less constrained 0.1 eV range would imply 
three orders of magnitude suppression. The asymmetry grows linearly with the mass of the 
heavy steriles. However, for masses larger than ~ 10 -100 GeV lepton number violating 
transitions via the Majorana mass could washout further the asymmetry, an effect that 
requires a refinement of the formulation to be taken into account. 
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3.3 Lepton asymmetries in the active sector 

The asymmetry generated ignoring the y evolution depends only on the Dirac-type phase, 
< 5 , appearing in W as we have seen. However when the evolution of the leptonic chemical 
potentials is included, other phases contribute to the total lepton asymmetry. We will 
perform a perturbative expansion to third order in the mixings of both V and W matrices. 
The result at finite t <C # 2 3 ( 0 2 3 ) 7 j” 1 can be written in the form: 

TfMW = J C pyl /cp (*) ( 3 - 45 ) 

Icp 

where all the four CP invariants appear, Iqp = j Jvk,- lj 2 \ } > given in eqs. ( 2 . 29 ). 

At finite t, the result for the functions Aj CP is well approximated by 

A (2 )(t) = ym(y 2 -yi )(l - 7vGi(t), 

1 v IN J 

Aj( 3 )(t) = -yiy 2 {yl - y\) (l - i 2 N G 2 {t), 

A V IN J 

A A 3 ) (t) = yiy 2 (l- 7ArG , 3 (t), 

A V 7 nJ 

Aj w (t) = 7i72 (l- Gn(t) - ^-G A 2 (t). ( 3 . 46 ) 

V 7 N J 27 AT 

where 77V = 7^ + 7^ and 77V = 27jv )) 37jv ^ w hQ e 

G\(t) = (e -72 * — e _7if ) Re [*J2o(Ai2, — A12, t) + 2 A„ J2oi(Ai2, — A12, t)] 

! 2 

+ 2 ^( — ^ ke lkt ^ e [' 4 2io(Ai2, — A12, t) (— 2 A„ + i{2^/k ~ 7 i — 72))] , 
z k =1 

( 3 . 47 ) 

and 

G 2 (t)=G 1 (% i=0 , ( 3 . 48 ) 

where we have defined 7* = y 2 7iv and A„ = v 2 — v\, and the result for G 3 (t), G4i(t), G42(t) 
are lengthier and reported in the appendix B. These results would get modified for 7,f 1 

had we included the non-linear terms that modify the rate of thermalisation at large times. 
In these equations there is an implicit expansion up to third order in 7j(uj)/A 1 / 3 when 
A 1 / 3 t lj while the terms 7 i(vi)t are resumed. 

In figure 3 we plot the functions A ( 2 )(t) and A 73) (t), which depend only on one 
neutrino mass difference. We show two physical situations: one with very degenerate 
neutrinos and the other with no strong degeneracies. 

These two invariants are the only ones relevant for the scenario that has been considered 
in most previous studies, where it has been assumed that only two sterile neutrinos have a 
role in generating the lepton asymmetry (see for instance [ 33 ] for a very recent analysis). 
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Figure 3: Functions A w(t) (left) and A (3 )(f) (right) assuming the rates are dominated by top 
-*1 L 1 

quark scattering, and taking y-ij \/2 = yi = 10 ~ 7 , for two choices of A Mf 2 = lGeV 2 (dashed) 
and AM 2 2 = 10 -6 GeV 2 (solid). t E w is the electroweak phase transition time, corresponding to 
T ew - 140 GeV. 


This is the situation in the limit of complete decoupling of IV 3 , ensured by the condition 
9i3 = 0, implying that only the invariants I ] ' and /]' ; survive. In [ 8 ] an approximate 
analytical solution was obtained, expanding in the yukawa’s, under the assumption that 
IAisI- 1 / 3 <C tEw “C 7/ 1 ■ In this limit, the result of eqs. (3.46) and (3.47) can be simplihed 
to 

Tr \n\{t E w) - — ((2/2 -y\)h -V2 1 ! }ym{y2 -y 1 ) U-7 n -™-• 

V 7 N) J- EW 

(3.49) 


Comparing with eq. (2.34), we see that the dependence on the yukawa’s is again that 
expected from a flavour invariant CP asymmetry. Using 


Im[J 20 (Ai 2 , -Ai 2 ,oo)] — —2 


/ 2\ 1/3 7 T 3 / 2 sign(Ai 2 ) 

V 3 y r[—1/6] |A 12 |2/3 ’ 


and 7 tv = | 7 at, and assuming the naive seesaw relations y\ = 2 ^ As v ° 2 ‘ Ml -, y\ 
we find: 


Tr[y](t EW ) ~ 10 


(3) v/m.mU / Ml \ 2 ' 3 
1 GeV 5 ' 3 IlAMfell 


(3.50) 

2 V&atmM 2 

i) 2 


(3.51) 


while for yf = y%/2 = 10 14 (that would correspond to light neutrino masses in the eV 
range and heavy ones in the GeV range) we would have 


Tr [y\{t E w) - 7 x 10 10 


4 2) - 21 { 3) 
|AM 2 2 (GeV 2 )| 2 / 3 ' 


(3.52) 


Even if the CP invariants are of 0(1), the asymmetry is too small unless there is a significant 
degeneracy between the two states [ 8 ]. It is important however to realise that the naive 
seesaw scaling is too naive and a full exploration of parameter space is necessary. 
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Figure 4: Functions -A r ( 3 ) (£) (left) and Aj w (t) (right) assuming the rates are dominated by 

top quark scattering, and taking yz/y/Z = Vi = 10 -7 , for three choices of [AM 2 2 ,AM 2 3 ] = 
[1, 2], [10 -6 , 2] and [10~ 6 , 2 x 10 -6 ] in GeV 2 (dashed, dotted and solid), tsw is the electroweak 
phase transition time, corresponding to Tew — 140GeV. 


In figure 4 we plot the functions A ( 3 ) (t) and A j w (t) . They depend on the two neutrino 
mass differences, so we show three examples here: one in which there are no degeneracies, 
one where there are two almost degenerate states, and the case where the three states are 
almost degenerate. As in the previous case we see a large enhancement when only one of 
the mass differences is small and a further enhancement when the two are small compared 
to the absolute scale. In the case of Aj w we find that there is a significant difference in 
the regime A-^'T -C 1 if we plot Aj w (t) truncated to the terms of O(yf). As we will see 
in the next section, the latter is much closer to the numerical result. The reason for this 
difference is that at small times, A l ^t <C 1, only some terms of order 0(yf ) are kept in 
eqs. (3.46), while there is a strong cancellation if all had been included. Note however that 
this effect is only important at times where the asymmetry is suppressed and seems to 
affect only Aj w . 

It is interesting to note that even though the dependence on the yukawas of the func¬ 
tions Aj CP (t) is different (fourth or sixth order), the maxima for all cases are roughly of the 

same order of magnitude. Note, however, that in the limit t 3 > 77 1 , only the contribution 

(3) 

of two invariants, Jw and , survive: 


lim Tr[/r](f) ~ 


AN 


lim 


27 tv 


Ap 3 3(t) 

Peq 


eq. (3.41) 


2 4 /3 vr 3/2 (3) ym (yl - y\) / 

3 1/3 r [-g] 1 |A 12 p/3 V 



7 N, 


(3.53) 


where we kept only the leading terms 0(y A ) proportional to if'' 1 and we have used the 
result of eq. (3.50). 

The first term in this expression corresponds to the expectation of [7], ie. the final 
asymmetry is proportional to that stored in the third sterile state, eq. (3.40), while the 
second term was missing in the simplified treatment of [7]. Note that they depend on 
different CP invariants. 
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Figure 5 : Left: full numerical solution (solid blue) and numerical solution neglecting non-linear 

terms (dotted green) for case 1, normalised to the invariant I\ , compared with the prediction, 

AA2)(t) (dashed red). Right: same for case 2 normalised to the invariant 7 -f 3 " 1 compared to A (3 )(f). 
■*1 ^1 
The parameters are the same as in figure 3 for the degenerate case. 

4. Numerical solution 

In order to check the accuracy of the analytical solutions presented in the previous sec¬ 
tion, we have solved the differential equations numerically. As shown in [28], the momen¬ 
tum dependence does not change significantly the results so we will consider the average- 
momentum approximation. 

In figures 5-6 we compare the analytical and numerical solutions for the functions 
Ai CP (t) in the highly degenerate case (the values of the mixing angles are of 0(1O~ 2 )) . In 
order to isolate the appropriate invariant we make the following choices: 

- _ ^2) 

• Case 1: 9^ = 9^ = 0 isolates I\ , 

• Case 2: 9i 3 = 9\2 = 0 isolates l[’\ 

• Case 3: 9\2 = 6 %3 = 0 isolates I^\ 

• Case 4: 9ij = 0 isolates Ij w . 

The numerical results normalised by the corresponding CP invariant are shown together 
with the predictions of the previous section. In the case of Jw-, we plotted the function 
Aj w keeping only the terms of 0(y 6 ) that is more accurate at small t and the full function 
at large t. The agreement in all cases is quite good. The differences observed at large t 
come from the non-linear terms in the equations. We also show the numerical results of 
the equations without them and find a very good agreement also at large t. Note that the 
approximation works well in the regime jt 3> 1 , that is in the strong washout regime of 
the fast modes. 

Numerically it is very hard to go to regimes where the ratios y/jAl 1 / 3 become very 
small, since the system becomes stiff. On the other hand, there is no reason why the 
perturbative solution is not accurate in such regime. We will therefore assume this to be 
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Figure 6: Left: full numerical solution (solid blue) and numerical solution neglecting non-linear 
terms (dotted green) for case 3, normalised to the invariant /^ , compared with the prediction, 
A ( 3 )(t) (dashed red). Right: same for case 4 normalised to the invariant Jw compared to Aj (t). 
The parameters are the same as in figure 4 for the double degenerate case. 


the case in the following section and use the perturbative solution to perform a scan of 
parameter space. 

5. Baryon asymmetry 

The observed baryon asymmetry is usually quoted in terms of the abundance, which is 
the number-density asymmetry of baryons normalised by the entropy. After Planck this 
quantity is known to per cent precision [34]: 

y“ p ~ 8.6(1) x 10~ n . (5.1) 

The lepton asymmetries in the left-handed (LH) leptons generated in the production of 
the sterile neutrinos are efficiently transferred via sphaleron processes [35] to the baryons. 
The baryon asymmetry is given by 

28 

y B = —Y B -l- (5.2) 

Since we have neglected spectator processes in the transport equations, the B — L asymme¬ 
try is related to the chemical potentials computed in the previous sections by the relation 

Y B -l = —^-TY[/x], (5.3) 

i r 4 #* 

where g * = 106.75 (which ignores the contribution to the entropy of the sterile states). 
Our estimate for the baryon asymmetry is therefore 

Ys ~ 3 x 10 -3 Tr[n(t)]\ tEw . (5.4) 

We have performed a first scan of the full parameter space of the model. Given the 
theoretical uncertainties mentioned in different sections of the paper, we have considered 
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as interesting the points that can explain the baryon asymmetry within a factor of 5. For 
this we have used the analytical solutions, even though in some regions of parameter space 
they will not be precise, since they are based on a perturbative expansion on the mixing 
angles of the matrices V and W. We have considered however a few cases where the angles 
are not small and we find that the analytical solutions differ from the numerical ones only 
in some global numerical factor of a few , but the time dependence is very similar. 

Even with an analytical expression the exploration of the large parameter space is a 
challenge. We have used the package Multinest [11, 12] to perform a scan on the Casas- 
Ibarra parameters [36], where the Yukawa matrix is written as 

_ ^2 

Y = -iUp MNS ^m hgh tR{zij) T \fM—. (5.5) 

miight is a diagonal matrix of the light neutrino masses and R is a complex orthogonal 
matrix that depends on three complex angles z t j. We fix the light neutrino masses and 
mixings to the present best fit points in the global analysis of neutrino oscillation data of 
ref. [37] and leave as free parameters: three complex angles, the three phases of the PMNS 
matrix, the lightest neutrino mass as well as the heavy Majorana masses that are allowed 
to vary in the range Mi G [0.1,100] GeV. In total thirteen free parameters. 

The scan searches for minima of the quantity | log 10 |Ys(tEw)/^ xp || (in the range 
< 1.5) and the MultiNest algorithm is optimised to sample properly when there are several 
maxima. For the determination of Yb we use the analytical results of the previous sections, 
for which the CP invariants are computed directly from the matrix elements of the V, W 
matrices that can be easily calculated by diagonalising the Yukawa mass matrix obtained 
in the Casas-Ibarra parametrization. Since the mechanism to work requires that at least 
one of the modes does not get to equilibrium before the electroweak phase transition we 
restrict the search to the range where one of the yukawa eigenvalues, 2 / 3 , is much smaller 
than the others and the following conditions are satisfied 

2/3 < 0.01Min[yi,y 2 ], ^ F 8 (|V)3| 2 + |Wi3 1 2 ) < 0.01I7 u (Tew)- (5 .6) 

*= 1,2 

Furthermore, since the kinetic equations neglect lepton number violating effects in the 
rates, we impose additionally the constraint 

/ m- \ 2 

* Ti<^H u {T EW ). (5.7) 

\J-EW / 

We first consider a case where one of the sterile neutrinos is effectively decoupled from 
baryon number generation, that we can assume to be IV 3 . This can be achieved with the 
choice of parameters: 


m 3 (i) = 0, z i3 = 0, R(zij) ->■ R(zij)(P), (5.8) 

for the IH(NH), where P is the 123 —> 312 permutation matrix (only necessary for the 
NH). With this choice, only the terms corresponding to the CP invariants I\" and /(' 
contribute. This case is the one that has been considered in most previous works on the 
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Figure 7: Points on the plane AM = M 2 — Mi versus Mi for which Yb > 1/5 x Y^ xp (blue), 
Yb > Fg Xp (green) and Yb > 5 x Y^ xp (red) for NH (left) and IH (right), with only two sterile 
neutrino species. 




subject [8, 9, 38, 28, 31], where the number of parameters is reduced to six: only one 
complex angle, two PMNS CP phases and two Majorana neutrino masses are relevant. 

It is believed that a large degeneracy of the two sterile neutrinos is needed to obtain the 
correct baryon asymmetry. In figure 7 we show the result of the scan under the conditions 
of eq. (5.8) on the plane AM 12 = M 2 — M\ versus Mi for normal and inverted orderings 
of the light neutrinos. The different colours correspond to values of Yb > 1/5,1, 5 x Y^ xp 
(blue,green,red). Successful leptogenesis is possible in a larger range of parameter space 
for IH than for NH. In the range shown our results agree reasonably well with those in 
ref. [39] for the IH, while the range for NH looks a bit smaller. We see that there are a 
significant number of points for which the degeneracy is mild for the IH. We have analysed 
more carefully some of these points by solving the full numerical equations. We find that 
even though these points correspond to cases where the angles in V, W are not small, the 
analytical and numerical solution agree very well and have the same t dependence as shown 
in figure 8. Note that the numerical solution is difficult at large times for non-degenerate 
solutions and the standard methods that we use fail. An optimised numerical method 
is needed to solve the stiffness problem and this will be studied elsewhere. It is very 
interesting to correlate the baryon asymmetry with observables that could be in principle 
measured such as the Dirac CP phase of the PMNS matrix, the amplitude of neutrinoless 
double beta decay or the active-sterile mixings that control the probability for the heavy 
sterile states to be observed in accelerators or in rare decays of heavy mesons. The effective 
mass entering the Ozx/5/3 decay is given by 

3 ' 3 ' J\yl 0l, PP (M) 

mpa = U 2 ei rm + Y, U 2 e{i+z) Mi MOv(lf} ' , (5.9) 

i= 1 i= 1 ' 

where JY[ 0u PP are the Nuclear Matrix Elements (NMEs) defined in [40] 3 . The first term 

3 The results for the NMEs computation in the interacting shell model [41, 42] are available in Appendix 
A of [40] 
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Figure 8: Comparison of the analytical (red-dashed) and numerical (blue-solid) solution for one 
of the points with mild degeneracy and Yb > Y^ xp , corresponding to log 10 (Mi(GeP)) = 0.9 and 
log 10 (AM (GeV)) = —0.92 and yukawa couplings yi = 1.3 x I0~ 6 ,y 2 = 9.8 x 10~ 9 . 


corresponds to the standard light neutrino contribution and the second is the contribution 
from the heavy states. U e i with i > 4 is the active-sterile neutrino mixing. 

In figure 9 we show the results for the active-sterile mixing as function of the sterile 
mass and compare them with present direct bounds and the prospects of SHiP[43] and 
LBNE near detector [44]. We show the result for Mi but the one for M 2 is almost identical. 
We see that most of the parameter space for successfull baryogenesis is not excluded by 
present constraints and that the active-sterile mixings tend to be larger for the IH. A 
sizeable region in the range of the GeV could be explored in the future experiment SHiP 
in the case of the IH and by LBNE near detectors. It is interesting to note that the less 
degenerate solutions can not have very small active-sterile mixing, as shown in figure 10 , 
where we plot the points on the plane €d eg = \M 2 — M\\/(M 2 + Ml) versus the active-sterile 
mixing in the electron flavour. The degeneracy can be lifted to some extent at the expense 
of larger yukawa couplings which also imply larger mixings. 

We have looked for direct correlations of the baryon asymmetry with the phases of the 
PMNS matrix. We have found that the distribution on the Dirac phase and the Majorana 
phase are flat. This is due to the fact that the complex angle can provide the necessary CP 
violation, even if the PMNS phases would vanish. The same is true for the effective mass of 
neutrinoless double beta decay, which depends on the Majorana phase. A dedicated scan 
is needed to quantify how the putative measurement of various observables could constrain 
the lepton asymmetry. This will be done elsewhere. 

In the general case, N ,3 is also relevant and the main difference with respect to the 
previous situation is that there is a significantly enlarged parameter space where degeneracy 
is not necessary. This was already found in refs. [46] for some points of parameter space. 
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Figure 9: Points on the plane |t/ e 4 1 2 (left), || 2 (middle), |[7^-41 2 (right) versus Mi for which Yb 
is in the range [1/5 — 1] x Y/] xp (blue) and [1 — 5] x Y/[ xp (green) for NH (up) and IH (down), with 
only two sterile neutrino species. The red bands are the present constraints [45], the solid black 
line shows the reach of the SHiP experiment [43] and the solid red line is the reach of LBNE near 
detector [44]. 



loflJU/ 


Figure 10: Points on the plane tdeg = M 2 2+M [ versus | t/e 4 1 2 for which the asymmetry is in the 
range [1/5,5] x Y^ xp in the range explored for IH. 


In figure 12 we show the points on the plane (AMi 2 ,Mi) for the general case. The active- 
sterile mixings are shown in figure 13. These mixings can be larger in this case, specially 
in the case of the NH. The SHiP prospects are therefore more promising in this context. 
As in the N = 2 case there is no direct connection between the asymmetry and the PMNS 
CP phases. On the other hand, the lightest neutrino mass is non-zero in this case, but the 
requirement that one yukawa needs to be significantly smaller than the others, eq. (5.6), 
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Figure 11: Distribution of mi for points that satisfy Yg > Y^ xp for the NH. 



Figure 12: Points on the plane AM = M 2 — Mi versus Mi for which Yg > 1/5 x Y^ xp (blue), 
Yb > Y/ xp (green) and Yg > 5 x Y/) xp (red) for NH (left) and IH (right), in the general case with 
three neutrinos. 


implies that the lightest neutrino mass must be small. In figure 11 we show the distribution 
of this quantity for those points that satisfy Yg > Y(/ xp in the case of NH (the IH being 
very similar). 

6. Conclusions 

We have studied the mechanism of leptogenesis in a low-scale seesaw model that is arguably 
the simplest extension of the Standard Model that can account for neutrino masses. For 
Majorana neutrino masses in the GeV range, sizeable lepton asymmetries can be generated 
in the production of these states some of which never reach thermal equilibrium before the 
electroweak phase transition. Lepton asymmetries are efficiently transferred to baryons 
via sphaleron processes. This mechanism was proposed in [7, 8] and studied in many 
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Figure 13: Points on the plane |C/ e 4 1 2 (left), \U^i\ 2 (middle), |C/,- 4 1 2 (right) versus M\ for which the 
Y b is in the range [1/5 — 1] x Y^ xp (blue) and [1 — 5] x Y/[ xp (green) for NH (up) and IH (down), 
with three sterile species. The red bands are the present constraints, the solid black line shows the 
reach of the SHiP experiment [43] and the solid red line is the reach of LBNE near detector [44] . 


works, but a full exploration of parameter space in the general case of three neutrinos 
is lacking. To this aim we have developed an accurate analytical approximation to the 
quantum kinetic equations which works both in the weak and strong washout regimes of 
the fast modes (there is always a slow mode that does not reach thermal equilibrium before 
the EW phase transition). It relies on a perturbative expansion in the mixing angles of the 
two unitary matrices that diagonalise the Yukawa matrix. This analytical approximation 
allows us to identify the relevant CP invariants, and explore with confidence the regime 
of non-degenerate neutrino masses which is very challenging from the numerical point of 
view. We have used this analytical solution to scan the full parameter space using the 
MultiNest package to identify the regions where the baryon asymmetry is within an order 
of magnitude of the experimental value. We have performed first a scan in the simpler 
setting where one of the sterile neutrino decouples, which reduces the parameter space, 
and is the approximation that has been considered in most previous works on the subject, 
for example in the so-called z^MSM. Although baryon asymmetries tend to be larger in the 
case of highly degenerate neutrinos, we find solutions with a very mild degeneracy that 
also correlate with a larger active-sterile mixing. These non-degenerate solutions appear 
for an inverted ordering of the light neutrinos. On the other hand we do not observe a 
direct correlation with other observables, such as the PMNS CP phases nor the neutrinoless 
double beta decay amplitude. 

We have also performed a scan in the full parameter space, with the only requirement 
that one of the yukawa matrix eigenvalues is very small, and that one mode will not reach 
equilibrium before the electroweak transition, for the washout not to be complete. The 
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main difference with the simpler case of two neutrinos is that the parameter space with 
successfull baryogenesis is significantly enlarged, in particular as regards non-degenerate 
spectra. Also the active-sterile mixings can reach larger values, particularly in the normal 
hierarchy case, improving the chances of future experiments such as SHiP or LBNE to find 
the GeV sterile neutrinos. There is much less difference in this case between normal and 
inverted neutrino orderings and also no direct correlation with the PMNS phases. On the 
other hand, the requirement of a small yukawa eigenvalue implies that the lightest neutrino 
mass cannot be large. 

A number of refinements are needed to improve the precision of the determination of 
the baryon asymmetry. First a more precise determination of the scattering rates of the 
sterile neutrinos is required. Most previous studies, and this one, have included only top- 
quark scatterings, but it has been pointed out recently that gauge scatterings are also very 
important. A correct treatment of these processes in the kinetic equations is necessary. 
Also the kinetic equations neglect effects of 0((Mi/T) 2 ). Such effects are not so small 
for masses in the GeV near the electroweak phase transition and their effect should be 
quantified. Finally, spectator processes and the asymmetries of fields other than the sterile 
neutrinos and LH leptons have not been taken into account in the kinetic equations. A 
proper treatment could easily bring corrections of 0(1). Finally, a more ambitious scan 
of parameter space should define more accurately the limits of eq. (5.6) for successfull 
baryogenesis. These effects will be studied in the future. 
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A. Results for the perturbative integrals 

A.l One dimensional integrals 

We just need them up to 0((3/a) 2 : 

(3 2 

J\(a\, P\,t) ~ Jio(ai,io) + P\Jn(ot\, to) + ^^12(«iTo) + tS.Ji(a\, /3i,t,to), 

(A.l) 
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with 


AJi(ai,/3i,Mo) = ^2 ^j<An(ai,Mo) - i ^ 


. ■“l ( n .ait 3 

/3? / e'V 1 e *^- 


+ O(r 4 ,t 0 - 4 ) 


n! 1 ccita n aR 2 


—n 


a*^ir fl +/ 3 i*o g*^-+/3it 


e 3 


= % 


aitl 


a\t 2 


(A-2) 


We can factor out the a dependence and define: 

J w {a,t) = |-t^(^Re[Jio(l,t|a| 1/3 )] + isign(a)Im[Ji 0 (l, f|a| 1/3 )]). (A.3) 

A.2 Two dimensional integrals 

We just need them up to 0((3/a): 

J2(&l:Pl,Ot2, /?2 ,t) — «/200(«li «2, to) + ^l«^21o(«l) a 2, to) + P 2 J 2 OI («1 1 «2 Ao) 

+ AJ 2 (ai, /3i, «2, /?2 , t, to), (A.4) 


where if a* / 0: 


AJ2(«1, Pi, CX- 2 -, 02 , t, to) — J\{oi 2 , P 2 , to) + * 


i° 2 > ~:hU,' 


e 3 




AJi(ai,Pi,t,t 0 ) 




.e 3 


— i it 


— i- 


a 2 J2i a ito a 2 J2i a it A 


(A.5) 


and for those terms where a* = 0 


e i°^ L +hto 


AJ 2 (« 1 , /?l, « 2 , ($ 2 , t, to) — \ J\{cy. 2 , 02, to) + i , 2 

a^o 

/ ft eEi^\ 


AJi(ai,/3i,t,t 0 ) 


«2 yjto y 

We can factorize the a-dependence: 

J 20 o(~Oi,a,t) = |^t_^Re[J 2 oo(—1, MH 1/3 )] + i sign(a)Im[j 20 o(-l, 1, t|a| 1/3 )]), 
J 20 i(-a,a,t) = yt_(Re[j 2 oi(-l, l,t|a| 1/3 )] + i sign(a)Im[j 20 i(-l, 1, f|a| 1/3 )]|, 
J 2 io(~ot,a,t) = |t-^Re[J 2 io(-l, l,t\a\ 1/3 )] + isign(a)Im[j 2 io(-l, 1, t|a| 1/3 )]), 


(A.6) 


(A.7) 


and reduce the integrals to the basic ones. 


- 29 



A.3 Three dimensional integrals 

We need the integrals up to 0(/3/a)° in this case. We can use the relation: 


^ ic 2X 3 


J 3 o(oii,oi 2 ,a 3 ,t) = Jio(ai,t)J 2 oo(oi 2 ,oc 3 ,t) — / dx e 3 J 1 o(a 1 , x)Jio(a3> x). 

Jo 


Since a* = 0 for the cases of interest, we can rewrite the result in terms of some basic 
integrals, I\ and I 2 : 

J 20 o(ai,a2,£) = t -^ r ^/i(|a 2 /ai|,sign(a 2 ),sign(ai),t|ai| 1/3 ) (A.9) 


rt 3 

ia 2 X 


dx e 3 Jio(ai, x)J\o ( 03 , x) = 


h (|«i/« 2 1, sign(ai), sign(a 3 ), sign(a 2 ), t\a 2 | 1/3 ) 


|aia 2 a3| 1,/3 


where 


h{r,si,s 2 ,t) = [ dx e lS2x3/3 J 10 (s 1 ,r 1/3 x) 
Jo 


(A.10) 


(A- 11 ) 


rt 

h(r,s 1 ,s 2 ,s 3 ,t) = / dx e lS3x3/3 J w (s 1 ,r 1/3 t)J 10 {s 2 ,(-s 3 /s 2 - s 1 /s 2 r) 1/3 t) (A.12) 
Jo 


B. Perturbative result for the invariants Jw and I 2 

The hnite t perturbative results proportional to the invariants Jw are given by the following 
expressions (we have used the property 7 \ oc 17 to simplify them): 

A Jw (t) = 7172 (l ~ G 41 (t) - ^G 42 (t), (B.l) 

V 7 N J 27 AT 

where 


_l)fc e -7 kt 


G 41 (t) = £(-!) 


Re CtijJ 2 o{^ij 1 ^ij 5 ^0) “ 1 “ bij J201 5 ^ij ? ^0) “ 1 “ C— (^) ^210 1 0) 


012 = * , 612 = 2 A.y, = -A„ + |( 27 fc - 71 - 72 ), 

ai3 = ~j bi 3 = 2vi, c[ 3 } = -vi - 2 ( 2 ^ - 71 ), 

023 = i, b 23 = -2u 2 , 4 ^ = ^2 + |( 27 fc - 72 ), 
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and A„ = v 2 - **i, A 7 = (72 - 71). 

G 4 2 ( t ) = Re d \ J3o(Ai2 — A13, —A12, A13, t ) + ^2^30(^12 — A13, A13, — A12, t ) 

+ ^3«^3o(Ai3, A12 — A13, — A12, t) + di J3o(Ai3, — A12, A12 — A13, t) 

+ d5«/3o(Ai2, — A12 + A13, — A13, t) + de J3o(Ai2, — A13, — A12 + A13, t) , 

(B.3) 

with 

d \ = Z \ + *y [ 2 A„ + *A 7 ] [ 2 n 2 - *(272 - 72)]e _72t , 

d 2 = ^2 + [ ~ 2*7 + *71] [2n 2 - *(272 - 72)] e _72t , 

d 3 = z 2 - *-y [ 2 v 2 + *72] [ 2 vi + *(27! -7 1 )]e _ ^ lt , 
d 4 = ^3 + *y [ 2 A„ - *A 7 ] [ 2 vi + *(271 - 7i)] e -71 *, 

d 5 = y (2n 2 + *72) e~' nt (2iA v + 271 - 71 - 72) - e _72t (2*A„ + 272 - 71 - 72) , 

d 6 = -y (2m - *71) e~' flt (2iA v + 271 - 71 - 72) - e _72 *(2*A„ + 272 - 71 - 72) • 

(B.4) 

(3) 

On the other hand, for the invariant ; 

Aj( 3) (i) = */i*/ 2 (l - =77) lNG 3 {t), (B.5) 

2 

G 3 (t) = ]T(-l) fc 

k = 1 

^ Re a' i jJ2o(Aij, —Aij,to) + 60(t)J2oi(Ajj, —A^-, to) + c^ k \t)J 2 io(Aij, — A^, to) 

i<j 

+Rc U>1 J3o(Ai2, — A12 + A13, —A 13 , t) + U>2</3o(Al2, — A 13 , — A 12 + A13, t) |, 

(B.6) 

with 

***i = 7 [2^271 + *7172], ***2 = ^ [-2^72 + *7172] , (B.7) 

and 

a' l2 = *72, &i 2 = 272^2 - **172 - **271, c/* ) = 572 ( - 2A,, + *(27^ - 72 - 71)), 

°13 = -*71, &13 = 2 7 i**i, c}^ = -571 ( 2 **i + *(27fc - 71)), 

«23 = *72, &23 = 272**2 , 4? ) = 272 (2n 2 + *(27fc - 72)) . 
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